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Abstract 

We describe the implementation of supersymmetric models in Pythia 8 including 
production and decay of superparticles and allowing for violation of flavour, CP, and R- 
parity. We also present a framework for importing generic new-physics matrix elements 
into Pythia 8, in a way suitable for use with automated tools. We emphasize that 
this possibility should not be viewed as the only way to implement new-physics models 
in Pythia 8, but merely as an additional possibility on top of the already existing 
ones. Finally we address parton showers in exotic colour topologies, in particular ones 
involving colour-epsilon tensors and colour sextets. 

1 Introduction 

The Standard Model (SM) of particle physics has been enormously successful in describing 
interactions between fundamental particles. The only experimentally unverified component 
of the SM is the Higgs boson which is thought to underlie the electroweak symmetry breaking 
(EWSB). However, experimental evidence for dark matter, neutrino masses and the theoret- 
ical requirement of naturalness appeals for a theory beyond standard model (BSM). Monte 
Carlo generators fulfil an important role, both in testing the SM to high precision and in 
testing the prediction of new theories, by providing a systematic procedure of comparing 
theoretical prediction to experimental observation. 

Pythia 8 [1] is a general-purpose Monte Carlo event generator [2] for a full simulation of 
high-energy collision events. It includes a comprehensive library of hard-scattering processes, 
particle decays, initial- and final-state parton-shower models [3lll], hadronization through 
string fragmentation [5] and models of beam remnants and multiple interactions [Sll^- It 
contains a native implementation of a wide variety of SM and BSM processes and also provides 
a standard interface [51[n] to external programs which may be used by a standalone generator. 

We describe here, the updates to the Pythia 8 event generator to include the popular 
BSM model of supersymmetry (SUSY), additions made to parton showers and hadronization 
algorithms to allow for exotic colour topologies and generic enhancements made to enable 
interfacing to parton-level BSM generators. In section [2l we describe the implementation of 
supersymmetric models in Pythia 8, including production and decay of superparticles. In 
section [3l we present a framework for importing generic new-physics matrix elements into 
Pythia 8, in a way suitable for use with automated tools. In section [H we discuss the 
treatment of parton showers in exotic colour topologies. Section [5] contains a brief summary 
and conclusions. 



2 Super symmetry in Pythia 8 



Supersymmetry (see |10j for a pedagogical introduction) is considered one of the best moti- 
vated extensions of the SM due to its abihty to address many outstanding theoretical and 
experimental issues. In particular, the Minimal Supersymmetric extension of the Standard 
Model (MSSM) is currently a popular candidate for a BSM theory. The MSSM extends the 
SM by the addition of one pair of SUSY generators which implies the presence of one super- 
partner to each SM state. The MSSM particle spectrum therefore has squarks (qi), sleptons 
(ii) and gauginos {B, and g) as the supersymmetric counterparts of quarks, leptons and 
gauge bosons respectively. The requirement of self-consistency of the theory via anomaly 
cancellation also demands two Higgs doublet fields Hu and H^. After electroweak symmetry 
breaking (EWSB), we are left with five Higgs degrees of freedom viz. the CP-even /iq and 
Hq, the CP-odd Aq and two charged Higgs bosons H^. The superpartners of the Higgses — 
the fermionic "Higgsinos" — mix with the gauginos to form neutralinos and charginos. In 
particular, the neutral Higgsinos {Hi and H2) mix with the neutral U{1) and SU{2) gauginos 
{B and W^) to form the mass eigenstates called the neutralinos = 1 — 4.) Similarly, 

the charged Higgsino mixes with the charged SU{2) gaugino to form charginos (Xo^ ~ 1)2.) 
The next-to-minimal supersymmetric extension of the SM (nMSSM) extends this scenario 
by adding one extra singlet Higgs field. This adds another member to the neutralinos and 
the neutralino mixing matrix is enlarged to 5 x 5. The current implementation of Pythia 8 
includes the nMSSM extension and allows processes with CP, flavour or R-parity violation. 

Pythia 8 uses the standard PDG codes for numbering the superpartners [TT] and the 
particle spectrum is read in via an SLHA file [12l[13]. We use the super-CKM basis (in the 
conventions of the SLHA2 [T3]) for describing the squark sector which allows non- minimal 
flavour violation. The mass-eigenstates of the squarks are then related to the left- and right- 
handed squarks via a 6 x 6 complex mixing matrix. Our implementation can therefore be 
used to study both CP violation and flavour violation in the squark sector. 
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The neutralino mixing matrix is a 4 x 4 (5 x 5 in the case of nMSSM) mixing ma- 
trix describing the transformation of the gauge eigenstate fermions {—iB, —iW^, Hi,H2) into 
the mass eigenstates (Xi) X2' Xa^ X4)- The two chargino mixing matrices U and V describe 
the diagonalization of the chargino mass matrix from the gauge eigenstates {—iW~^,H~^) to 
ixtiXt)- Supplementary conventions for vertices and most of the cross-section formulae are 
taken from [H], as detailed below. 

2.1 Couplings 

Pythia 8 reads particle masses and mixing matrices via the SUSY Les Houches Accord 
(SLHA2) framework [13]. (Read-in of SLHAl spectra [12j is also supported, but mixing the 
two standards is strongly discouraged, as the internal translation from SLHAl to SLHA2 
has only been designed with the original SLHAl in mind.) The raw data read in by the 
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SusyLesHouches class is accessed by the CoupSUSY class which uses the information to con- 
struct all the SUSY couplings. The couplings are defined according to [2] for all cases except 
for couplings of superparticles to Higgs bosons which are defined according to |15j . 

The running of electroweak and strong couplings is carried over from the corresponding 
one-loop calculations in the Standard Model. The GAUGE block can be used to set the boundary 
values of all three SM couplings at the SUSY breaking scale. By default, the masses of W 
and Z are assumed to be the pole masses and are used to calculate the on-shell value of 
sin^ 6\y = 1 — m^r/m^. If externally provided in the SLHA file, the value of sinOw can be 
set to the running value using the flag SUSY: sin2thetaWMode = 2 (see the Pythia 8 HTML 
user reference included with the code [1]). The ratio of the two Higgs vacuum expectation 
values (tan /3) is read in from the low scale value provided by the MINPAR and EXTPAR blocks. 
The default value of the Higgs mixing angle (an) is set to the SM limit (/3 — it/2) which is 
then overwritten by the contents of the HMIX block. 

Since the SLHA interface has been extended and can now be used to pass information on 
any new particles and decays [9] , the presence of the MODSEL block is used as an indicator of 
SUSY models and Pythia 8 will initialize the CoupSUSY class only if this block is present. 
Skipping the MODSEL block is acceptable for Les Houches Event files (LHEF) as long as the 
user supplies an external decay table for all required cascade decays. 

2.2 R-parity violation 

The most general MSSM superpotential allows both lepton and baryon-number violating 
processes. This is generally avoided by demanding invariance under an R-parity defined as 
(|_1)3B-L+2S'_ ppom this definition, all SM particles are even whereas all superpartners are 
odd under R-parity. A well known consequence of this is that the Lightest SUSY particle 
(LSP) must be stable. A neutral, weakly interacting LSP can therefore be a good candidate for 
dark matter. However, the imposition of R-parity can be considered an aesthetic requirement 
rather than a consistency requirement and possible R-parity violating interactions, if present, 
can be probed by collider experiments. We therefore include R-parity violating production 
and decay processes in our implementation. 

In SLHA conventions, the R-parity violating superpotential is given by 

Wrpv = fiiH^Li + ^X,jkL,L,Ei, + KjkLiQjDj, + ^X^^^iDpl (2) 

The /i-type terms correspond to bi-linear R-parity violation which causes a mixing between 
the leptons and neutralinos/charginos. The A and A'-type terms lead to lepton number vio- 
lation whereas A"-type terms lead to baryon-number violation. The current implementation 
does not include the effects of the bi-linear term. The R-parity violating couplings Xijk are 
antisymmetric under i -H- j. Therefore only couplings for i > j are read and the rest are 
set by the symmetry property. Similarly, X'-jj, is antisymmetric under j ■v^ k and hence 
only couplings with j > k need to be provided. This implementation includes in particular, 
the resonant production of a squark via A"-type couplings which can be probed at hadron 
collider experiments. The changes made to showering and hadronization to account for the 
non-standard colour structure from such terms will be explicitly described in section HI 

2.3 Cross Sections 

The current implementation of SUSY includes all leading-order (LO) 2 — )• 2 production pro- 
cesses with gluinos, squarks, charginos, and neutralinos in the final state and also 2 — t- 1 — t- 2 
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Subprocess class 


processname 


Chargino and neutralino production 


qqbar2chi0chi0, 

qqbar2chi+-chi0, 

qqbar2chi+chi-. 


Gaugino squark production 


qg2chi0squark, 
qg2chi+-squark. 


Gluino production 


gg2gluinogluino, 
qqbar2gluinogluino. 


Squark-gluino production 


qg2squarkgluino 


Squark-pair production 


gg2squarkantisquark, 

qqbar2squarkantisquark 

qq2squarksquark 


RPV resonant squark production 


qq2antisquark 



Table 1: List of SUSY production processes. In all cases, charge conjugate processes are 
turned on by default. 



resonant production of squarks via baryon number violating couplings. All available SUSY 
processes can be turned on using SUSY : all = on. Individual subprocesses can then be se- 
lected based on the final state by setting SUSY: idA = PDGcode and SUSY: idB = PDGcode. If 
only idA is provided, all processes with that particle in the final state are turned on. Alterna- 
tively, one or more production processes can be turned on using the string SUSY: processname 
= on, again with SUSY : idA and SUSY : idB providing a further level of subprocess selection. 
The available subprocess classes are listed in Table [TJ 

The squark-antisquark and squark-squark production processes include contributions from 
EW diagrams and their interferences. To estimate the size of these contributions, and/or for 
purposes of comparison to other codes that do not include them, the cross sections can be 
restricted to include only the strong-interaction contributions, using the following flags: 

• qqbar2squarkantisquark: onlyQCD = true. 

• qq2squarksquark: onlyQCD = true. 

The baryon number violating coupling X'^-^ if present, can induce resonant squark produc- 
tion via the process djd^^ — >• u* which produces a resonant up-type antisquark or via Uidj — )• 
or Uidk — d*j which produce a down-type antisquark. The expression for an up-type squark 
production process is 

* jk i' 

The expression for down-type squarks is similar, taking into account the symmetry prop- 
erty X'-jf, = —^ikj- We implement this production process as qq2antisquark and the charge 
conjugate process {qiQj — >• ijk) is included by default. 

The supersymmetric Higgs sector is identical in many ways to the Two-Higgs Doublet 
Model. The Higgs production processes have already been implemented in Pythia 8 in 
the SigmaHiggs class. The production of the Higgs bosons can be accessed by including 
the switch HiggsBSM: all=on. For specific Higgs processes, please refer to the HTML user 
reference included with the code [1]. 
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2.4 Sparticle Decays 



SUSY Particle decays are handled by the class SUSYResonanceWidths. The user can choose 
to read in decay tables via SLHA or use the decay widths calculated by Pythia. As a default, 
Pythia does not calculate the decay width if a table is externally supplied. Note, however, 
that while Pythia's internal treatment can include sophistications such as matrix-element- 
based phase-space weighting and running widths, channels read in from an SLHA decay table 
will be decayed purely according to phase space, with no matrix-element weighting. The 
internal treatment should therefore be preferable, in most cases, and an option for overriding 
the automatic read-in of decay tables is provided, by setting the flag SLHA :useDecayT able = 
false, see sec. 13.11 

The decay of a particular particle may be turned off manually using the standard Pythia 8 
structure PZJGcode :mayDecay = false or by setting its width to zero in the SLHA decay 
table. In the former case, the particle will still be distributed according to a Breit-Wigner 
distribution with non-zero width, whereas it will always be assigned its pole mass in the latter. 

Individual decay modes may be switched on/off using the standard Pythia 8 methods, 
documented in the section on "The Particle Data Scheme" in the program's HTML documen- 
tation [1]. We discuss ways to switch modes on/off using SLHA decays tables in section [3?T1 

The internal treatment of 2-body decays is so far restricted to on-shell particles. A mecha- 
nism for effectively generating 3-body decays via sequences of 1 — )■ 2 decays involving off-shell 
particles is foreseen as an update in the near future (and will be announced in the Pythia 8 
update notes). An equivalent mechanism is already implemented in Pythia 8, e.g., for 
h —7- ZZ decays for light Higgs bosons. 

Currently the following R-parity conserving two-body decays are implemented: 

• g ^ m-j 

• x° ^ kij, xp, x]w- 

• xt ^ mj, kh, xjz, x°w+ 

• qi^Qjxt ijxt, qjZ, qjW+ 

Besides these, we also include two-body R-parity violating decays of squarks via A' {q — t- 
Iq') and A"-type couplings {q — t- q'q"). We also include the three-body decays of neutralinos 
through A"- type couplings via an intermediate squark |16j . For certain final states in three 
body decays, partial decay via sequential two-body decays may also be kinematically allowed. 
In this case, we demand that only the off-shell components of the matrix element-squared are 
allowed to contribute to the three-body decay width. Any interferences between the off-shell 
and on-shell components are also turned off. The two-body sequential decays then proceed 
as normal. 

The Higgs boson running widths are calculated in the associated classes ResonanceH for 
CP even {ho, Hq) and the CP odd (^o) Higgses, and ResonanceHchg for charged Higgses(i?^). 
By default, the Higgs decay tables are not overwritten even if they are read via SLHA because 
Pythia 8 performs a more accurate phase space calculation than the flat weighting that is 
performed for decay widths read in via SLHA. The decays of Higgses into SUSY particles will 
be included in a future update. 
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3 Interfacing Generic BSM Models 



The simplest way of implementing a new model may often be to just inherit from SM or 
BSM processes that have already been implemented in Pythia 8, modifying and generalizing 
them as appropriate, as described in the section on "Semi-internal Processes" in the main 
Pythia 8 documentation [1]. 

Alternatively, Pythia 8 can read in parton-level events generated by external matrix- 
element event generators |17ft21j . using the Les Houches Event File (LHEF) format [8ll22j. 
If the events contain new particles, so-called QNUMBERS blocks described in sec. 13.11 

below, can be used to add information on the quantum numbers of new particles, and SLHA 
decay tables [12] may also be provided. (For SUSY models, in addition, complete SLHA 
spectra can be given, as discussed in sec. [2j) The encoding of colour flow is then particu- 
larly important, for the events to be showered and hadronized correctly. Some pedagogical 
examples, with illustrations, are given in the original LHA paper (22], and further explicit 
examples with colour-epsilon and colour-sextet structures are given in section [J] below. The 
LHEF paper [8] describes how to encode this in an LHE file, with examples of correct LHE 
files available, e.g., in Pythia's examples/ directory. 

When reading events from LHE files, the BSM/SLHA information may either be en- 
closed within the LHE file (preferred), or provided in a separate file. In the former case, the 
BSM/SLHA information should be included in the <header> part of the LHE file inside 
an <slha> tag [9]. In the latter case, a separate BSM/SLHA file may be specified using the 
Pythia 8 command SLHA : file = fileName. The mode SLHA : readFrom gives the user some 
additional control over whether and from where BSM/SLHA information is read in. It should 
normally be left at its default setting, but can optionally be used either to switch off SLHA 
read-in entirely, or to force read-in from a specific file: 

SLHA: readFrom = # do not read BSM/SLHA information at all 

= 1 # (default) read in from the <slha> . . . </slha> block of a LHEF, 

# if such a file is read during initialization, and else from 

# SLHA: file 

= 2 # read in from SLHA: file 

The framework described in section 13.21 represents a third option which combines fea- 
tures from both of the two possibilities above. It allows parameters and matrix-element 
code to be imported directly from external packages, to generate semi- internal processes in 
Pythia 8 (i.e., without an intermediate LHE file) in a fully automated and generic way, as 
long as the final-state parton multiplicity does not exceed the limitations of Pythia's in- 
ternal hard-process phase-space generator [1]. A working interface between Pythia 8 and 
MadGraph 5 |20j has been constructed along these lines, for 2 — )■ 2 processes, and will be 
reported on separately. Here, we focus on the Pythia 8 side of the interface. 

The interface basically consists of two parts: 1) information about particles and couplings 
using a generalized SLHA format (section l3.ip . and 2) accessing that information from within 
a semi-internal Pythia 8 process (section 13. 2p . 

3.1 Information about New Particles 

Information about particle quantum numbers, masses, couplings, and decays, can be given 
in an ASCII file, using a generalization of the SLHA [12] and BSM-LHEF [9] formats, whose 
name is provided to Pythia by setting the word SLHA: file = fileName. 
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3.1.1 QNUMBERS 



The SLHA file should contain a QNUMBERS block [9] for each state not already associated with 
an ID code (a.k.a. PDG code, see [Ulllj for a list) in Pythia 8. For a hypothetical electrically 
neutral colour-octet self-conjugate fermion (a.k.a. a gluino) that we wish to assign the code 
7654321 and the name "balleron", the structure of this block should be 

BLOCK QNUMBERS 7654321 # balleron 

1 # 3 times electric charge 

2 2 # number of spin states (2S+1) 

3 8 # colour rep (1: singlet, 3: triplet, 6: sextet, 8: octet) 

4 # Particle/Antiparticle distinction (O=own anti) 

For a non-selfconjugate particle, separate names can be given for the particle and its 
antiparticle. For a heavy up-type quark, 

BLOCK QNUMBERS 8765432 # yup yupbar 

1 2 # 3 times electric charge 

2 2 # number of spin states (2S+1) 

3 3 # colour rep (1: singlet, 3: triplet, 6: sextet, 8: octet) 

4 1 # Particle/Antiparticle distinction (O=own anti) 

Note that the name(s) given after the # mark in the block definition are optional and 
entirely up to the user. If present, they will be used, e.g., when printing out event records 
with Pythia's event. list () method. 

The SM quantum numbers given in the QNUMBERS blocks are required by Pythia 8 for 
QED and QCD showering, and for colour-flow tracing. (Currently, Pythia does not make 
use of the spin information.) As a rule, we advise to avoid clashes with existing ID codes, 
to the extent possible in the implementation. A useful rule of thumb is to only assign codes 
higher than 3 million to new states, though one should be careful not to choose numbers 
larger than a 32-bit computer integer can contain, which puts a cap at ~ 2 billion. 

3.1.2 MASS 

The file should also contain the SLHA block MASS, which must, as a minimum, contain one 
entry for each new state, in the form 

BLOCK MASS 

# ID code pole mass in GeV 

7654321 800.0 # m(balleron) 
8765432 600.0 # m(yup) 

In principle, the block can also contain entries for SM particles. Here, some caution and 
common sense must be applied, however. Allowing SLHA spectra to change hadron and/or 
light-quark masses in Pythia 8 is strongly discouraged, as these parameters are used by the 
parton-shower and hadronization models. Changing the 6-quark mass, for instance, should 
ideally be accompanied by a retuning of the b fragmentation parameters. Since this is not the 
sort of question a BSM phenomenology study would normally address, by default, therefore, 
Pythia 8 tries to protect against unintentional overwriting of the SM sector via the fiag 
SLHA:keepSM, which is on by default. To be more specific, this fiag causes particle data 
(including decay tables, see below) for ID codes in the ranges 1-24 and 81 - 999,999 to be 
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ignored. Notably this includes (ID 23), (ID 24), and t (ID 6). The SM Higgs (ID 25), 
however, may still be modified by the SLHA input, as may other particles with ID codes in 
the range 25 - 80 and beyond 1,000,000. If you switch off this flag then also SM particles are 
modified by SLHA input. 

Alternatively, the parameter SLHA:minMassSM, with default value 100.0 GeV, can be spec- 
ified to allow any particle with ID code below 1,000,000 to be modified, if its default mass in 
PYTHiAlies below some threshold value, given by this parameter. The default value of 100.0 
allows SLHA input to modify the top quark, but not, e.g., the and bosons. 

3.1.3 DECAY 

The file may also include one or more SLHA decay tables New BSM particles without 
decay tables will be treated as stable by Pythia 8. For coloured states, this may result in 
errors at the hadronization stage, and/or in the possibly unintentional production of so-called 
i?-hadrons [25] , with a reasonably generic model for the latter available in Pythia 8 [26] . On 
the other hand, a redefinition of Pythia's treatment of the decays of SM particles, like Z^ and 
may be undesirable, since Pythia's internal treatment is normally more sophisticated 
(discussed briefly in sec. 12. 4p . Thus, again, caution and common sense is advised when 
processing (B)SM particles through Pythia, with the protection parameters SLHA:keepSM 
and SLHA : minMassSM described above also active for decay tables. An option for overriding the 
automatic read-in of decay tables is also provided, by setting the flag SLHA :useDecayT able 
= false. 

The format for decay tables is [12] 

# ID WIDTH in GeV 

DECAY 7654321 2 . 034369169E+00 # balleron decays 

# BR NDA IDl ID2 ID3 
9.900000000E-01 3 6 5 3 # BR( -> t b s ) 
l.OOOOOOOOOE-02 3 4 5 3 # BR( -> c b s ) 

Note that the branching ratios (BRs) must sum up to unity, hence zeroing individual 
BRs is not a good way of switching modes off. Instead, Pythia 8 is equipped to interpret a 
negative BR as a mode which is desired switched off for the present run, but which should be 
treated as having the corresponding positive BR for purposes of normalization. 

Finally, a note of warning on double counting. This may occur if a particle can decay via 
an intermediate on-shell resonance. An example is — )• qiq2Q3QA which may proceed via 
WW followed by W — > qq'. If branching ratios for both —?■ WW and — )• Qi 92(73^4 
are included, each with their full partial width, a double counting of the on-shell — >■ WW 
contribution would result. (This would also show up as branching ratios summing to a value 
greater than unity.) Such cases should be dealt with consistently, e.g., by subtracting off the 
on-shell contributions from the — t- qiq2qzq4: partial width. 

3.2 Accessing the Information from a Semi-Internal Process 

Already the original SLHAl [12] allowed for the possibility to create user-deflned blocks, be- 
yond those deflned by the accord itself. The only requirement is obviously that the block 
names already deflned in the accord(s) should not be usurped. The SLHA interface in 
Pythia 8 will store the contents of all blocks, both standard and user-defined ones, as internal 
vectors of strings. 
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By default, Pythia's internal BSM implementation only extracts numerical content from 
those blocks it recognizes (i.e., the standard SLHA 1&2 blocks and QNUMBERS), and uses those 
to initialize its couplings and particle data arrays. However, generic methods are also provided, 
that can be used access to the contents of any block, whether standard or user-defined, from 
inside any class inheriting from Pythia's SigmaProcess class (i.e., in particular, from any 
semi-internal process written by a user), through its SLHA pointer, slhaPtr, by using the 
following methods: 

bool slhaPtr->getEntry (string blockName, doublefe val) ; 

bool slhaPtr->getEntry (string blockName, int indx, doublefe val); 

bool slhaPtr->getEntry (string blockName, int indx, int jndx, doublefe val); 

bool slhaPtr->getEntry (string blockName, int indx, int jndx, int kndx, doublefe val); 

This particular example assumes that the user wants to read the entries (without in- 
dex, indexed, matrix-indexed, or 3-tensor-indexed, respectively) in the user-defined block 
blockName, and that the entry value, val, should be interpreted as a double. In fact, the 
last argument is templated, and hence if anything other than a double is desired to be read, 
the user has only to give the last argument a different type. Since the user presumably knows 
what type of content his/her own user-defined blocks contain, this solution allows the content 
to be accessed in the correct format, without Pythia needing to know what that format is 
beforehand. If anything goes wrong (i.e., the block does not exist, or it does not have an 
entry with that index, or that entry cannot be read as a double), the method returns false; 
true otherwise. This effectively allows input of completely arbitrary parameters using the 
SLHA machinery, with the user having full control over names and conventions. Of course, it 
is then also the user's responsibility to ensure complete consistency between the names and 
conventions used in the SLHA input, and those assumed in any user- written semi-internal 
process code. 

Note also that the special SLHA block SMINPUTS (containing SM parameters [12j) will al- 
ways be accessible through the methods above, regardless of whether a corresponding SLHA 
block has been read in or not. The SMINPUTS block is initialized starting from PYTHIA's own 
internal default values, with subsequent modifications as dictated by updates to PYTHIA's 
particle and parameter databases before initialization and/or by SLHA read-in. This func- 
tionality is intended to give a generic BSM implementation access to the SM parameters 
contained in SMINPUTS in a universal way. 

To give a specific example, the interface to Madgraph 5 was structured in the following 
way. Among the possible output formats available for matrix elements in Madgraph 5, one 
is a mode called pythiaS. When invoked, this mode writes out the corresponding matrix 
element(s) in exactly the format required by Pythia 8's semi-internal process machinery. 
The resulting code can therefore be imported directly into Pythia 8, and Madgraph even 
provides explicit instructions and a Makefile for doing precisely that. In general, however, such 
matrix elements may contain parameters that refer, e.g., to couplings in a model unknown 
to Pythia. A central question was therefore how to provide information on such parameters 
at runtime, in a sufficiently generic manner. The solution is that Madgraph writes out the 
relevant parameters as custom-made SLHA- like blocks in a BSM/SLHA file included together 
with the matrix-element code. It then also inserts appropriate calls to slhaPtr->getEiitry () 
in the cross-section expressions, so that each parameter can be retrieved when needed, without 
any user intervention required. 

Note that this entirely circumvents a particularly troublesome issue that before was 
thought to make any truly universal "BSM Accord" impractical, the problem of agreeing 
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on a common standard for names and parameters for completely arbitrary models. In the 
Madgraph-Pythia interface, it is sufficient that Madgraph itself assigns some unique 
names and contents to each block. It has complete freedom in choosing which conventions to 
use, as long as it correctly inserts the corresponding readEntryO calls in its matrix-element 
output. This effectively generalizes the SLHA data structure to apply to completely general 
BSM models. 

The interface has been tested by the authors (in collaboration with our Madgraph col- 
leagues) to work for importing both a few trivial examples of models, such as a 4-generation 
model, to more exotic ones, such as a model with colour-sextet diquarks and one with a 
baryon-number-violating vertex. Showers in such topologies are the topic of section HI 

4 Showers and Hadronization in Exotic Colour Topologies 

In this section, we describe Pythia's treatment of QCD radiation in topologies containing 
colour-epsilon tensors (section 14. ip and ones involving particles with colour-sextet quantum 
numbers (section l4.2p . This applies regardless of whether the event is generated as an internal, 
semi-internal, or LHEF process. We also comment briefly on hadronization aspects, pointing 
out relevant sources of further information. 

As an aid to implementations using LHEF, a few examples of how to arrange Les Houches 
colour tags in colour-epsilon and colour-sextet cases are given in Figs. [T] and [TcJ respectively. 
For completeness, we also show the status and ID codes, and the mother information, for each 
particle. These are explained in more detail in [22], where also more illustrations (including 
both standard and baryon-number violating ones) can be found. 

For completeness. Fig. [id] shows a situation which PythiaIs not yet capable of handling. 
The illustration shows a complicated baryon-number-violating cascade decay of a hypothetical 
fourth-generation top quark (assigned ID code 8) involving both supersymmetric, fourth- 
generation, and SM particles, to produce a situation with a total of three colour-connected 
baryon-number-violating vertices. At the moment, Pythia's junction fragmentation model 
|27j is at most capable of handling up to two connected colour junctions (specifically, single 
junctions and junction-antij unction systems), hence only if a g ^ qq branching in the shower 
happens to break up the triple-junction system into smaller colour-singlet subsystems would 
Pythia's fragmentation model be able to deal with it. 

A somewhat less pathological case in which multi-junction topologies may result is if a 
single baryon-number violating vertex becomes colour-connected to both of the junctions in 
the (baryon) beam remnants. This may happen some small fraction of the time through 
multiple parton interactions. In such cases, the following error message will be printed and 
the generation of the event restarted, 

Error in ColConf ig: : insert : junction topology too complicated; too many junction legs 

4.1 Colour-Epsilon Topologies 

For colour topologies involving the epsilon tensor in colour space (i.e., colour topologies with 
non-zero baryon number) we first consider the example oit^qq in the RPV-SUSY model. 
The Lagrangian for the UDD-type interaction terms is 

C = -A^^-.e'"- {uUd^)fPRdl + d^A,{u^)\PRdl + dMu^)\PRdl + h.c) (4) 
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Figure 1: Illustration of the assignment of Les Houches colour tags in four different exotic 
colour topologies. Lines corresponding to colour (anticolour) tags are drawn above (below) the 
propagators. Top Left: s6 — )• t* — )• txi- Top Right: e~e+ — ;> Z° — ;> (t* -)• sb){i sb). Bottom 
Left: production of a colour-sextet particle, assigned the fictitious ID code 6000001; the 
negative anti-colour tag (drawn below the sextet propagator) is interpreted as an additional 
(positive) colour tag. Bottom Right: A complicated baryon-number-violating cascade decay 
(of a hypothetical fourth-generation fermion) producing a total of three colour-connected 
baryon-number-violating vertices; such topologies (with three or more interconnected colour 
junctions) cannot yet be handled by Pythia's string fragmentation model [27j . 
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Figure 2: Gluon emission from RPV vertices with e-tensor. 



To extract the behaviour of the radiation function, we look at the ratio of exact matrix 
element for ti?(pi) — )■ d{p2)s{p-i) + g{q) via A3^2 to the matrix element for — )■ d{p2)s{p^) 
and retaining only the parts that are soft- or collinear-singular (i.e., which diverge for one or 
more q ■ Pi — )• 0). Since momentum is explicitly conserved in the shower branching process, 
the pre- and post-emission momenta must be related by 



Pl=P2+P7i=P2+Pi + q 

h p\ = m\ = s the invariant mass of the decaying squarl 
The Born-level matrix element squared is given by: 

1 2 _ 1 2 



(5) 



\M< 



01 



\>:k2V{Nc-i)\s 



(6) 



Three diagrams (shown in Fig [2]) contribute to the process where one gluon is emitted 
from this configuration. The matrix element corresponding to this process i.e. t/?(pi) — )• 
d{p2)s{pz)g{q) is denoted by Mi and, for massless decay products (p^ = p\ = 0); is given by 



I Ml 



2gl\\"\\Nc 



1 



Nc 
1 

+ - 



ly.Cps 
(pi ■P2) 

1 \{pi- q){p2 ■ q) 



{P2 ■ q) ^ {P3 ■ q) ^ 



+ 



{Pi ■ Ps) 



{pi ■ q){P3 

X 



+ 



{P2 ■ Ps) 



q) {p2-q){p3-q) 



(P3 • q) {P2 ■ q) Nc 



1 



+ Y 



where 



X 
Y 
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Pi ■ q 

{P2 ■ P3){'m^ 



{P2 ■ Ps) 

pi-q) 



{P2 ■ q) 



+ 



{P3 ■ q) 



(pi ■ q) 



{pi-q){p3-q) {pi-q){p2-q) {p2-q){p3-q) 



(pi ■ qf 



(7) 



(8) 



The antenna pattern represented by equation ([7]) can be characterized as follows: the terms 
on the second line represent three soft-eikonal dipole factors (see, e.g., [2]), one for each of 
the three possible two-particle combinations. The factor l/{Nc — 1) in front of the dipole 
factors implies that the normalization of each of these eikonals is half as large as that of the 
eikonal term in an ordinary qq antenna, see, e.g., |28H30j . The two first terms on the last line 
of equation ^ correspond to additional purely collinear singularities for each of the quarks. 
The factor l/(A''c — 1) is here absent; the collinear singularities have the same strength as those 
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of an ordinary qq antenna. The X and Y terms contain subleading-color and a quasi-collinear 
term for the decaying t, respectively. 

The eikonal terms (including the leading part of the Y term, oc m? /{pi ■ q)"^) agree with 
the expression in |16p31j. which is used to generate radiation for this type of colour topology 
in Herwig [32]. (Note that, in the Herwig implementation, the pattern is generated using 
ordinary full-strength radiation functions, by selecting randomly between each two-particle 
combination, thereby reproducing the full pattern when summing over events |16p31j.) 

For the implementation in Pythia, we have chosen a different strategy, as follows. First, 
using momentum conservation, we may rewrite the antenna pattern above to only contain 
the final-state particle momenta, 

VI-P2 ^ Pi -Ps ^ P2 -Pz ^ 2 



{pi-q){P2-q) {pi-q){p3-q) iP2-q){p3-q) pi-q 

This reduces the eikonal part of expression to a single antenna between the two final-state 
quarks, plus subleading leftover terms. The eikonal and the collinear terms then correspond 
exactly to the standard radiation pattern from a qq dipole with an extra term of 0{^). 

For the present work, we therefore take the radiation pattern of a standard-strength dipole 
spanned between the two final-state quarks as our starting point. Using Sij = 2pi -pj, this 
radiation function is given by [28] 

= 8.a.cJ^ + 4^ + 4^) . (10) 

\Mz^qq\^ \S2qS3q SS3q SS2q J 



The Pythia showers are not based directly on equation (jlOp . but rather on Altarelli-Parisi 
(AP) splitting kernels, which partition the radiation pattern onto two terms, each of which is 
governed by the q ^ qg splitting function, 

Pq^q.iz) = Cf^-^ , (11) 

with z the energy fraction retained by the quark after emitting the gluon. The energy fractions 
of the final-state quarks, 2 and 3, are defined as in [33], 

X2 mj-S3q X3 rn\ - S2q 

Z2 = = ; Z3 = = — 5 . (12) 

X2 + Xq m{ + S2q X3 + Xq mf + S^q 

The expression actually used in the Pythia showering is the sum of the AP contributions, 

1*^'!' 8™,fQ2> + ^) , (13) 



I^Op V S2q S3q 

The full matrix-element ratio, |Mip/|Mop, as well as the various approximate forms 
discussed here, are illustrated in Fig. [21 with the mass of the decaying t arbitrarily set to 
mi = 300 GeV. On the left-hand pane, we show the size of the radiation function (without 
the overall factor SnasCp) as a function of the opening angle between the final-state gluon 
and one of the quarks, for a fixed (~soft) gluon energy Eg = 10 GeV. On the right-hand 
pane, we show the dependence on energy, for a fixed (~collinear) opening angle 6qg = 20°. 
The bottom row shows the ratio of each approximation to the matrix-element result. 

The thick solid (blue) line represents the full t* qq + g matrix element, equation ([7]). 
For comparison, the thin solid (red) line shows the pattern obtained for a standard dipole. 
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Figure 3: Illustration of radiation functions for gluon emission in q* — t- qiq2 decays for soft 
(left) and coUinear (right) gluons. 



equation ()10p . The dashed (brown) curve shows the Pythia approximation to the dipole 
pattern, given by the sum of the AP splitting kernels in equation (|13p . Finally, the light 
dot-dashed (green) curve shows the eikonal approximation to the matrix element, used by 
Herwig. 

In the soft limit (left-hand pane of Fig. [3|), all the expressions agree in the two extremal 
points, in which the gluon is both soft and collinear. For wide-angle soft emissions, e.g. at 
90° opening angle, the standard dipole pattern (as well as its DGLAP variant) overestimate 
the full matrix element by up to a factor ~ 1.5. That is, the Pythia shower will generate 
slightly too many soft wide-angle gluons. By contrast, as would be expected in the soft limit, 
the eikonal approximation works well for all opening angles. 

In the collinear limit (right-hand pane of Fig. [3|), the x axis is now the gluon energy, with 
the opening angle held fixed. All the expressions again agree for small gluon energies, in the 
double soft- and collinear limit. For intermediate gluon energies, the standard dipole pattern 
(as well as its DGLAP variant) again slightly overestimate the full matrix element, while they 
again agree with the matrix element in the hard collinear limit, on the right-hand edge of the 
plot. The eikonal, however, does not include the collinear-singular terms on the last line of 
equation ([7|) and hence does not reproduce the rise of the other curves in the hard-collinear 
limit. 

In summary, our shower model will slightly overestimate the total amount of radiation, in 
particular at large angles, while the Herwig model underestimates hard-collinear radiation. 
We therefore regard the two as complementary. We note that the neglected terms could 
still subsequently be incorporated into Pythia 8 as a matrix-element correction |341 l35j. 
presumably mostly relevant if S-violating processes should indeed be observed in nature. 

For the case x? ~^ QQQ^ the corresponding expression is similar to equation ([7]) with three 
half-strength eikonals between the quark from the neutralino-quark-squark vertex and the 
two quarks from the RPV vertex [16]. Besides these, each quark has the corresponding full 
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collinear singularity. Herwig treats this situation by randomly connecting each quark in 
the final state to either of the other two quarks. We have chosen instead to implement it 
in Pythia 8 as three genuinely half-strength dipoles spanned between the three final-state 
quarks. 

For the case of three-body RPV gluino decay, g — )• qqq, only the resonant parts, g — )• 
qq* — )• qqq, have so far been implemented in Pythia 8, cf. section 12.41 For future off-shell 
contributions, the emission structure of the non-resonant parts will be obtained from the 
relative strengths of the intermediate off-shell g — t- qq* contributions. 

In all cases, the subsequent hadronization phase makes use of Pythia's ability to handle 
string topologies including colour junctions [27], and hence issues such as baryon-number flow 
should be treated at least semi-realistically, allowing studies at the individual-particle level. 

4.2 Colour-Sextet Particles 

Within the leading- A'^c dipole approach to radiation adopted in Pythia, we represent a colour- 
sextet charge as the sum of two colour-triplet charges, in much the same way as octet charges 
(e.g., gluons) are represented as the sum of a triplet and an antitriplet charge. Each triplet 
charge is independently colour-connected to an antitriplet charge. Hence a sextet may be 
colour-connected either by a "double bond" to an anti-sextet (in an overall singlet 66 config- 
uration), or by two "single bonds" to two independent antitriplet charges, depending on the 
colour fiow in the event. Each such "bond" is interpreted as an ordinary QCD dipole, with 
the sextet end treated as a massive quark. 

At the technical level, we note that the Les Houches colour-tag standard was not originally 
designed to deal with sextet colour configurations. This is easy to remedy, however. Since 
a sextet never carries an anticolour, its anticolour tag is effectively available for use. To 
distinguish an additional colour (i.e., a sextet) from the ordinary anticolour (octet) case, we 
adopt the convention that a negative anticolour tag is interpreted as an additional colour, and 
vice versa for anti-sextets, as was illustrated in Fig. [Tel This appears to violate no present 
use of the standard (negative colour tags were so far never used in practise, as far as we are 
aware) . 

We note that a more complete treatment of the radiation and phenomenology of sextet 
diquarks was published while this manuscript was in preparation, see [36]. 

5 Summary and Conclusions 

We describe the implementation of Supersymmetry in the Monte Carlo event generator 
Pythia 8. We use the generic super-CKM basis of |13] which allows CP and fiavour vio- 
lation in the squark sector. We also allow R-parity violation in production processes and 
decays and the extension of the MSSM to the nMSSM. The current implementation includes 
all pair-production processes with gluinos, squarks, neutralinos or charginos in the final state. 
We also implement the resonant production of squarks via R-parity violating vertices that 
can be relevant at a hadron collider like the LHC. Two-body decays of all SUSY particles 
(except the Higgs sector) and R-parity violating decays of neutralinos via the A" couplings 
have been implemented. The Higgs decays will be implemented as a part of a future update. 

We also describe the enhancements made to the SLHA interface to allow external programs 
to pass non-standard information blocks to Pythia 8. The modifications provide a mechanism 
for so-called semi-internal processes in Pythia 8 to access all information read in via the 
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SLHA interface. This interface can therefore be used for implementation of generic BSM 
models without requiring a previous agreement on standardization of names and parameters. 

Finally, we have commented on how Pythia 8 handles showering in non-standard colour 
topologies, such as the epsilon topologies encountered in R-parity violating models and in 
sextet di-quark ones. 
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Appendix A: Test Sparticle Spectrum 



All validations have been performed using point SPSla (mSUGRA parameters mo = 250, 
^^1/2 = 100,^0 = 0,/i > and tan/3 = 10). However, since the masses and mixings of super- 
particles at low scale depend on renormalization group running, we give here the complete list 
of masses and mixing matrices used in our validations. The following spectrum was generated 
using SoftSUSY 2.0.5 [37] 
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Appendix B: Validated cross sections 



We present validated cross sections for point SPSla. All sparticle decays are turned off. 
The non-default parameters used were chosen mostly for simplicity, and to enable direct 
comparison with both the Pythia 6 |38p39j and XsuSY [T3] implementations: 

PDF:pSet = 8 (CTEQ6L1) 

SigmaProcess:factorscale2 = 4 (v^) 

SigmaProcess:renormScale2 = 4 (v^) 

SigmaProcess:alphaS value = 0.1265 

SigmaProcess:alphaSorder = 1 
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